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We study the nature of the phase transition in the fuhy frustrated simple cubic 



i-Q I lattice with the XY spin model. This system is the Villain's model generalized in 



three dimensions. The ground state is very particular with a 12-fold degeneracy. 

Previous studies have shown unusual critical properties. With the powerful Wang- 

c/3 , Landau flat-histogram Monte Carlo method, we carry out in this work intensive 

C^ ' simulations with very large lattice sizes. We show that the phase transition is clearly 

I ■ of first order, putting an end to the uncertainty which has lasted for more than twenty 

o 
o 



years. 
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I. INTRODUCTION 

One of the most fascinating tasks of statistical physics is the study of the phase transition 
in systems of interacting particles. Much progress has been made in this field since 50 
years. Finite-size theory, renormalization group analysis, numerical simulations, ... have 
contributed to the advance of the knowledge on the phase transition. Exact methods have 
been devised to solve with mathematical elegance many problems in two dimensions. But as 
improvements are progressing, new and more complicated challenges also come in from new 
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discoveries of materials and new applications. Renormalization group which has predicted 
with success critical behaviors of ferromagnets has many difficulties to deal with frustrated 
systems. Numerical simulations which did not need huge memory and long calculations 
for simple systems require now new devices, new algorithms to improve convergence in 
systems with extremely long relaxation time, or in systems whose microscopic states are 
difficult to access. One class of these systems is called 'frustrated systems' introduced in the 
seventies in the context of spin glasses. These frustrated systems are very unstable due to the 
competition between different kinds of interaction. However they are periodically defined (no 
disorder) and therefore subject to exact treatments. This is the case of several models in two 
dimension^i, but in three dimensions frustrated systems are far from being understood even 
on basic properties such as the order of the phase transition (first or second order, critical 
exponents, ...). Let us recall the definition of a frustrated system. When a spin cannot fully 
satisfy energetically all the interactions with its neighbors, it is "frustrated". This occurs 
when the interactions are in competition with each other or when the lattice geometry 
does not allow to satisfy all interaction bonds simultaneously. A well-known example is the 
stacked triangular antiferromagnet with interaction between nearest-neighbors. 

The frustration in spin systems causes many unusual properties such as large ground 
state (GS) degeneracy, successive phase transitions with complicated nature, partially dis- 
ordered phase, reentrance and disorder lines. Frustrated systems still constitute at present 
a challenge for investigation methods. For recent reviews, the reader is referred to Ref. y. 

In this work, we are interested in the nature of the phase transition of the classical XY spin 
model in the fully frustrated simple cubic lattice (FFSCL) shown in Fig. [TJ Although phase 
transition in strongly frustrated systems has been a subject of intensive investigations in 
the last 20 years, many aspects are still not understood at present. One of the most studied 
systems is the stacked triangular antiferromagnet (STA) with Ising^, XY and Heisenberg 
spins^i^. The cases of XY (A^ = 2) and Heisenberg (A^ = 3) STA have been intensively 
studied since 1987"^— , but only recently that the 20-year controversy comes to an end.— ""— 
For details, see for example the review by Delamotte et al^. 

The paper is organized as follows. Section II is devoted to the description of the model 
and the technical details of the Wang-Landau (WL) methods as applied in the present paper. 
Section HI shows our results. Concluding remarks are given in section IV. 




FIG. 1: Fully frustrated simple cubic lattice. Discontinued (continued) lines denote antiferromag- 
netic (ferromagnetic) bonds. 

II. MODEL AND WANG-LANDAU ALGORITHM 



We consider the FFSCL shown in Fig. [TJ The spins are the classical XY model of 
magnitude S = 1. The Hamiltonian is given by 



(id) 



(1) 



where Sj is the XY spin at the lattice site i, J2(i,j) is made over the NN spin pairs Sj and Sj 
with interaction Jij. Hereafter we suppose that Jij = —J {J > 0) for antiferromagnetic (AF) 
bonds indicated by discontinued lines in Fig. [1], and Jij = J for ferromagnetic (F) bonds. 
This model is a generalization in three dimensions (3D) of the 2D Villain's model^^ which 
has been extensively studied^^"— : every face of the cube is frustrated because we know that 
a plaquette is frustrated when there is an odd number of AF bonds on its contour-i^. To 
describe the model, let us look first at the xy plane (Fig. [T]). There, all interactions are F, 
except one AF line out of every two in the y direction. The same is for the yz [zx) plane: 
one AF line of every two in the z (x) direction. Note that there is no intersection between 
AF lines. Each plane is thus a 2D Villain's model. 

Let us recall some results on the present model. For the classical XY model on the FFSCL, 
the GS is 12-fold degeneracy with non collinear spin configurations^^. For convenience, let 
us define the local field acting on the spin Sj from its neighbors as h,j = J2j Jij^j- The GS 
can be calculated by noticing that the local field is the same at every site and is equal to 
Ihl = 2^3 so that^ 



2(S2 + S3 + S4) 



(2) 



he = 2(Si + S3 - S4) (3) 

h7 = 2(Si - S2 + S4) (4) 

hg = 2(Si + S2 - S3) (5) 

where the factor 2 resuhs from the symmetric neighbors lying outside the cube and J^- = ±1 
depending on the bond has been used. Putting into square these equahties and using 
h^ = 12, S^ = l(i = 1,...,8), one has three independent relations which determine the 
relative orientation of every spin pair 

S2 • S3 + S3 ■ S4 + S2 ■ S4 = (6) 

-Si ■ S3 + S3 ■ S4 + Si ■ S4 = (7) 

Si ■ S2 + S2 ■ S4 - Si ■ S4 = (8) 

There are 12 solutions of these equations which can be described as follows2^>22. Consider 
just one of them shown in the upper figure of Fig. H) On a yz face, the spins (displayed 
by continued vectors) on a diagonal are perpendicular, i.e. Si±S2, SrXSs. In addition, the 
orthogonal dihedron (Si, S2) makes an angle a = arccos( ^^^ ) with the dihedron (S7,S8). 
On the other yz face the spins (displayed by discontinued vectors in Fig. |2]) are arranged 
in the same manner: Ss-LSg, S4_LS3 and the dihedron (S4, S3) makes an angle a with the 
dihedron (S5, Sg). Note that the dihedron (S5, Sg) makes a turn angle /3 = 7r/4 with respect 
to the dihedron (S7,S8), and that the sum of the algebraic angles between spins on each 
face of the cube is zero. 

There is another choice shown in the lower figure of Fig. [2] where every thing is the same 
as described above except /3 = — 37r/4. One has therefore two configurations with the choice 
of the yz faces. If one applies the same rule for the spins on the xy faces and the xz faces, 
one obtains in all 6 configurations. Finally, together with their 6 mirror images, the total 
degeneracy is 12.— 

The above description of the GS shows an unusual degeneracy which can help to under- 
stand the first-order transition shown below by relating this system to a Potts model where 
the GS degeneracy plays a determinant role in the nature of phase transition. Note however 
that all simulations have been carried out with the initial Hamiltonian ([T]), no GS rigidity 
at finite temperature (T) has been imposed. 
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FIG. 2: Two of the 12 ground-state configurations of the fuhy frustrated simple cubic lattice. The 
numbers indicate the spins at the sites defined in Fig. 1. The angle a is a = arccos(^-'^). Upper: 
/3 = 7r/4, lower: /3 = -37r/4. 

The first investigation of the nature of the phase transition of this model by the use of 
Metropohs Monte Carlo (MC) simulations has shown a second order transition with unusual 
critical exponents^^. However, MC technique and computer capacity at that time did not 
allow to conclude the matter with certainty. Recently, Wang and Landau^S proposed a MC 
algorithm for classical statistical models which allowed to study systems with difficultly 
accessed microscopic states. In particular, it permits to detect with efficiency weak first- 
order transitions^--'^ The algorithm uses a random walk in energy space in order to obtain 
an accurate estimate for the density of states g{E) which is defined as the number of spin 
configurations for any given E. This method is based on the fact that a flat energy histogram 
H{E) is produced if the probability for the transition to a state of energy E is proportional 
to g{Er\ 

We summarize how this algorithm is implied here. At the beginning of the simulation, 
the density of states (DOS) is set equal to one for all possible energies, g{E) = 1. We begin 
a random walk in energy space (E) by choosing a site randomly and flipping its spin with a 
probability proportional to the inverse of the temporary density of states. In general, if E 
and E' are the energies before and after a spin is flipped, the transition probability from E 



to E' is 

p{E^E') = mm[g{E)/g{E'),l]. (9) 

Each time an energy level E is visited, the DOS is modified by a modification factor / > 
whether the spin fiipped or not, i.e. g{E) — )■ g{E)f. At the beginning of the random walk, 
the modification factor / can be as large as e^ ~ 2.7182818. A histogram H{E) records the 
number of times a state of energy E is visited. Each time the energy histogram satisfies a 
certain "flatness" criterion, / is reduced according to / — )■ ^/f and H{E) is reset to zero 
for all energies. The reduction process of the modification factor / is repeated several times 
until a final value /final which close enough to one. The histogram is considered as fiat if 

H{E)>x%.{H{E)) (10) 

for all energies, where x% is chosen between 70% and 95% and {H{E)) is the average 
histogram. 

The thermodynamic quantities^i^ can be evaluated by {E"') = 

^j:EE''9iE)exp{-E/kBT), C, = ^^^^, (M") = ^Ee M^g{E) expi-E/ksT), and 
X = - — ^ y , where Z is the partition function defined hy Z = J2e Qi^!) exp{—E/kBT). 
The canonical distribution at any temperature can be calculated simply by 
P{E,T) = y{E)exp{-E/kBT). 

In this work, we consider a energy range of interest^>^ (i^min, -E'max)- We divide this 
energy range to R subintervals, the minimum energy of each subinterval is -E^in f*^^ ^ — 
1,2,...,R, and maximum of the subinterval i is -En^ax ~ -^min + 2A£', where AE can be 
chosen large enough for a smooth boundary between two subintervals. The WL algorithm 
is used to calculate the relative DOS of each subinterval (-Emin; -^max) with the modification 
factor /final = exp(10^^) and fiatness criterion x% = 95%. We reject the suggested spin fiip 
and do not update g{E) and the energy histogram H{E) of the current energy level E if the 
spin-fiip trial would result in an energy outside the energy segment. The DOS of the whole 
range is obtained by joining the DOS of each subinterval (-Emm + ^E, -E^ax ~ AE). 

III. RESULTS 

We used the system size oi N x N x N where A^ varies from 24 up to 48. We stop at 
A^ = 48 because, as seen below, the transition at this size shows a definite answer to the 



problem studied here. Periodic boundary conditions are used in the three directions. J = 1 
is taken as the unit of energy in the following. 

We show in Fig. |3]the magnetization and the susceptibility and in Fig. IHthe energy per 
spin and the specific heat, for A^ = 24. All these quantities show a transition with an second- 
order aspect. However, we know that many systems show a first-order transition only at 
very large sizes. This is indeed the case. The energy histograms obtained by WL technique 
for three representative sizes A^ = 24, 36 and 48 are shown in Fig. El As seen, for A^ = 24, 
the energy histogram, though unusually broad, shows a single peak indicating a continuous 
energy at the transition as observed before in Fig. |H The double-peak histogram starts 
already at A^ = 36 and the dip between the two maxima becomes deeper with increasing 
size, as observed at A^ = 48. We note that the distance between the two peaks, i. e. 
the latent heat, increases with increasing size and reaches ~ 0.03 for A^ = 48. This is 
rather large compared with the value ~ 0.009 for A^ = 120 in the XY STA— "— i^ and with 
0.0025 for A^ = 150 in the Heisenberg case^i. We give here the values of Tc for a few sizes: 
Tc = 0.68080 ± 0.00010, 0.67967 ± 0.00010 and 0.67919 ± 0.00010 for A^ = 24, 36 and 48, 
respectively. 

Note that the double-peak structure is a sufficient condition, not a necessary condition 
in old- fashion MC simulations (i. e. not WL method), for a ffist-order transition. This 
is because in old-fashion MC simulations performed at a given T, we often encounter the 
situation where, at the transition, the system is spatially composed of two (or more) parts: 
the ordered phase with energy Ei and the disordered phase with energy E2. Since in old- 
fashion MC simulations, we make histogram from the total system energy, i. e. Ei + E2, 
the histogram will record the 'average' energy Ei + E2, therefore the double peak structure 
will not be seen. Such a coexistence at any time of the ordered and disordered phases in 
some ffist-order transitions makes it impossible in old-fashion MC simulations to detect two 
peaks. In our present WL flat-histogram method, the double-peak structure is obtained 
from the DOS histogram which gets rid of the problem of spatial coexistence of the two 
phases discussed above. Therefore, the double-peak structure is a necessary condition for a 
ffist-order transition as it should be. 

Let us say a few words on the correlation length. It is known that the correlation length is 
finite at the transition point in a ffist-order transition. For very strong ffist-order transitions, 
this correlation length is short so that the ffist-order character is detected in simulations 
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FIG. 3: Magnetization and susceptibility versus T for N =24. 
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FIG. 4: Energy per spin and specific heat versus T for N =24. 



already at small lattice sizes. For weak first-order transitions, the correlation length is very 
long. To detect it one should study very large lattice sizes as in the present paper. Direct 
calculation of the correlation length is not numerically easy. Fortunately, one has other 
means such as the WL method to detect more easily weak first-order transitions. 

To close this section, let us emphasize two points: i) First, the first-order transition 
observed here may come from the fact that the GS of the present XY FFSCL model has a 
12-fold degeneracy which is reminiscent of the 12-state Potts model. In three dimensions, 
the latter model has a first order transition. Note however that this conclusion is not 
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FIG. 5: Energy histogram for several sizes A^ = 24, 36, 48 at Tc indicated on the figure. 

always obvious because the continuous degrees of the order parameter mask in many cases 
the symmetry argument based on discrete models^l, ii) Second, some other XY frustrated 
systems such as the FCG^^, HCP— and hehmagnetio^^ antiferromagnets show also a first- 
order transition in MC simulations. Though the nonperturbative renormalization group has 
been extensively used for the STA case- due to its long-lasting controversy, we believe that 
the other cases are worth to study in order to verify that the validity of that theory is not 
limited to the STA but is universal for frustrated systems of vector spins. 

IV. CONCLUDING REMARKS 



Using the powerful WL fiat histogram technique, we have studied the phase transition 
in the XY fully frustrated simple cubic lattice. The technique is very efficient because it 
helps to overcome extremely long transition time between energy valleys in systems with a 
first-order phase transition. We found that the transition is clearly of first-order at large 
lattice sizes in contradiction of early studies using standard MC algorithm and much smaller 
sizes^^. The result presented here will serve as a testing ground for theoretical methods such 
as the renormalization group which still has much difficulty to deal with frustrated systems^. 
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